Quasi-stationary states in the self-gravitating sheet model 
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We study quasi-stationary states (QSS) resulting from violent relaxation in the one-dimensional 
self-gravitating "sheet model" , revisiting in particular the question of the adequacy of the theory of 
Lynden-Bell (LB) to describe them. For "waterbag" initial conditions characterized by a single phase 
space density, the prediction of this theory is, in this model, a function of only one parameter, which 
can conveniently be chosen to be the ratio of the energy to that in the degenerate limit. Studying a 
class of such initial conditions in which the shape of the initial waterbag is varied, we find that the 
LB predictions are reasonably good always in the low energy region, while at higher energies (i.e. in 
the non-degenerate limit) they are generally not even qualitatively correct, although certain initial 
conditions can still be found where they are as good as at low energy. We find notably that, in line 
with what has been observed by Levin et al. in some other models, when LB theory does not work 
the QSS are always characterized by the presence of a degenerate core, which these authors explain 
as the result of dynamical resonances. In short LB theory appears to be a good approximation 
only when violent relaxation is sufficiently "gentle" , and otherwise a degenerate core-halo structure 
results. 

PACS numbers: 



I. INTRODUCTION 

The rich statistical mechanics of long-range interacting 
systems has been a subject of active study in recent years 
(for a recent review see e.g. 0]). As for self-gravitating 
systems, such systems have been understood to give rise 
generically to non-equilibrium states which evolve only 
on time-scales which diverge with the number of parti- 
cles. The degree to which such "quasi-stationary" states 
(QSS) can be understood, and their properties predicted, 
by a statistical approach is a question which is inevitably 
posed. In this context a theory originally formulated by 
Lynden-Bell in the astrophysical context in the sixties 
and which has been applied also in the study of two 
dimensional vortices has seen revived interest in re- 
cent years. Study notably of a one-dimensional (ID) toy 
model, the Hamiltonian Mean Field (HMF) Model which 
describes particles on a ring interacting by a cosine poten- 
tial, showed that this theory can predict sometimes very 
accurately the properties of these states (see, e.g., |^ 
and references therein), and more generally manages to 
capture the qualitative dependence of the QSS on initial 
conditions. While it is clear that the "LB theory" is not 
entirely adequate in general, these studies suggest that 
the basic physical principle behind it — maximization 
of an entropy subjected to the constraints appropriate 
to Vlasov dynamics — is, at the very least, a reference 
point for understanding the out of equilibrium dynam- 
ics of these systems. This contrasts strongly with the 
view of this theory in the (original) context of the astro- 
physical literature, where it has simply been discarded 
as a completely inadequate, and basically irrelevant, the- 
ory 1^ . One recent study |Q of three dimensional (3D) 
self-gravitating systems concludes, however, that LB the- 



ory may indeed be relevant also to this case. This study 
shows that in a certain limited range of initial condi- 
tions the LB theory predicts well the density profiles of 
QSS, and proposes an alternative theory to explain their 
properties in the regime where LB no longer works well. 
The same authors have shown that the same statements 
apply both to plasma systems ^ and two-dimensional 
self-gravitating system [|^, and, in a very recent article 
PI , have used the alternative theory to account for QSS 
in the HMF model. 

In this paper we study these issues in the so-called self 
gravitating "sheet model" (SGS) of particles in one di- 
mension attracted by forces independent of separation. 
Our main goal is to characterize more precisely than has 
been done previously the degree of validity of the LB 
theory in this model, which is one of the canonical toy 
models for the study of such systems, and to determine 
whether the properties of the QSS can be characterized 
in a simple manner and perhaps understood when the LB 
theory does not apply. That this theory does not provide 
an adequate theory of QSS in the SGS model is clear from 
the earliest studies of this issue p^-[T^ which indeed used 
this model to probe the possible validity of LB theory for 
3D gravitating systems. More recently a study of these 
questions in the SGS model has been reported by Yam- 
aguchi who finds reasonable agreement with LB in 
a certain range of initial conditions, and, like in the work 
of Levin et al. mentioned above, proposes a modification 
of it to account for the QSS observed in other cases. We 
will compare in detail our results to these previous works. 

Studies of the SGS model in the astrophysical con- 
text go back at least as far as that of [|6|, and there 
have been numerous studies of it also in the statisti- 
cal mechanics literature in the decades since. Many of 
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these studies focussed on the question of relaxation to 
the thermal equilibrium of the model, for which the ex- 
act expression was first derived by Rybicki . That this 
relaxation, like in other long-range systems, takes place 
very slowly, on a timescale which diverges with the num- 
ber of particles, has been clear since the earliest works, 
but the precise time scale and parametric dependences 
thereof have been the subject of considerable study and 
even some controversy (see e.g. ]l8|-^ and references 
therein). In a recent work pT[ on this question, we have 
established clearly that the relaxation time from a range 
of initial conditions depends linearly on the number of 
particles N ^, while also showing a strong dependence 
on the intermediate QSS state (or states). Besides the 
early and more recent studies cited above which consider 
the QSS attained on the shorter mean-field time scales 
(i.e. through violent relaxation) and LB theory, there 
are also studies ^ which argue that the assumption 
that QSS always result from mean-field dynamics may 
not be always correct: starting from certain initial con- 
ditions the initial phase of relaxation is observed to lead 
to phase space densities which have large holes which ro- 
tate in phase space, which persist on the time scales of the 
simulations. In our analysis below we will examine this 
question carefully, as it is clearly of central importance 
to understand whether the formation of a QSS is indeed 
a good description of the outcome of violent relaxation if 
one is comparing with a theory which, by construction, 
assumes such a outcome. 



The article is organized as follows. First we will start in 
section m with the definition of the model, and its numer- 
ical integration. In the following section III we review the 



theory of violent relaxation of Lynden-Bell and describe 
our calculations of the predictions for the density profiles, 
and velocity and energy distributions . We will also in- 
troduce a simple set of "order parameters" which we use 
to characterize the QSS. In section [V we describe the 
specific class of initial conditions which we investigate. 
In section ^ we report our numerical results, comparing 
them to the theoretical LB predictions. In the following 
section we confront our results with two proposals which 
have been made in the recent literature to explain the 
properties of QSS when the LB is clearly inadequate. We 
also discuss our results briefly in the light of the kinetic 
theory for coUisionless relaxation developed in [ p6[ (and 
references therein) . In our conclusions we summarize our 
findings and conclusions, and suggest some directions for 
further investigation in both ID and 3D self-gravitating 
systems. 



II. THE SELF-GRAVITATING SHEET MODEL 

We consider identical particles of mass m in one dimen- 
sion which are mutually attracted by a force independent 
of their separation, i.e., the force on a particle i due to a 
particle j is 
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-gTO^sgn(xi - Xj 



where 5 is a couping constant. If the particles in one 
dimension are considered as infinitely thin parallel sheets 
in three dimensions interating by 3D Newtonian gravity, 
it is simple to show that gm = 27rI]G where G is Newton's 
constant and I] is the mass per unit surface area of the 
sheets. In a system of a finite number of such particles 
the total force acting on the i*^ particle at any time may 
be expressed as 



(1) 



where denotes the number of particles on the right 
of i*^ particle and Nt_ for the left. 

The fact that the force is thus constant other than 
when particles cross leads to one of the very nice features 
of this toy model; its numerical integration requires only 
the solution of algebraic (quadratic) equations to deter- 
mine the time of the next particle crossing. This means 
that the only limit on the precision of integration is that 
of the machine in solving such equations, and that no 
numerical parameters need be introduced. Another sim- 
plification comes from the fact that, in one dimension, 
the crossing of two particles without discontinuity in the 
velocities is, up to labelling of the particles, equivalent to 
an elastic collision in which particles exchange velocities. 
If we are not interested in following the trajectories of 
individual particles, we can thus consider the system as 
consisting of particles on which the forces are constant 
in time [and given by the initial value of (|l|)], and which 
undergo elastic collisions when they collide. The opti- 
mal way to treat this kind of problem is, as has been 
pointed out and discussed in detail in [^7|, by using a 
so-called "heap-based" algorithm, which uses an object 
called a "heap" to store in an ordered way the next cross- 
ing times of the pairs. This algorithm requires a num- 
ber of operations of order log(A^) to determine which of 
the — 1 pairs crosses next. Given that the number of 
crossings per particle per unit time grows in proportion 
to N, the simulation time thus grows in proportion to 
N'^ \og{N). As is common practice we will use the total 
energy (which is conserved in the continuum model) as a 
control parameter. For the longest simulations we report 
the error in total energy of the order of 10~^%. 



^ This result is consistent, notably, with the analysis of |22| which 
argues that a timescale linear in A'^ arises because of "resonances" 
present in spatially inhomogencous QSS in ID systems, but not 
in spatially homogeneous QSS which occur in ID systems such as 
the HMF model, where a faster scaling with N is indeed observed 
(see, e.g., I^sj). 



HI. PREDICTIONS OF LYNDEN-BELL 
THEORY 

In this section we very briefly recall the basics of the 
theory of Lynden-Bell, and describe how we calculate 
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its predictions for different quantities, in the case of wa- 
terbag initial conditions. 

A statistical theory to describe the stationary states 
arising from violent relaxation through mean-field forces 
has been proposed by Lyndcn-BcU in 1967 Such 
states were proposed to arise from the relaxation of the 
coarse-grained phase space density to that derived by 
maximizing the entropy derived for the latter by "count- 
ing" the (fine-grained or "microscopic" ) phase space con- 
figurations states consistent with the conservation laws 
imposed by the collisionless (Vlasov) dynamics. For the 
case of an initial "water-bag" phase space density, i.e., in 
which the microscopic phase space density has the same 
value everywhere it is non-zero, these conservation laws 
simply require the conservation of the phase space vol- 
ume "occupied" by this constant density, /o say. The 
calculation of the entropy is then equivalent to that for 
identical particles with a "fermionic" exclusion, and gives 
(in one dimension) 



nlnn -|- (1 — n)\n{l — n) >dxdv (2) 



where n = f / fa, and / is the coarse-grained phase space 
distribution in the macrocell at {x,v). Maximization of 
(|) gives 



f{x,v) 



fo 



I _(_ g/3(£(a:,-u)-M) 



(3) 



where e{x,v) = + 'fiix) denotes the energy density of 
phase-space element at (x, v). The constants /3 and /i are 
Lagrange multipliers associated with the conservation of 
the total mass M, and total energy E, of the system: 



M = 



f{x, v)dxdv 

,2 



E = I lC- + ^)f{x,v)dxdv., 



(4) 
(5) 



where (^(a;) is the mean field potential generated by the 
mass density p(x). Except in the degenerate and non- 
degenerate limits, corresponding to /3 — >■ c» and /3 — >■ 
respectively, it is not possible to solve these equations 
analytically to derive (/3,/i) for any given M, E and /q. 
It is, however, straightforward to do so numerically, as 
described in Appendix ^ (see also Q ) . 

We note that, although the prediction of LB theory for 
a water-bag initial condition depends in general on the 
three parameters M, E and /o, for the SGS model there is 
only one additional dimensionful quantity relevant in the 
continuum limit, the coupling g. Thus units can always 
be chosen so that two of M, E and /o are fixed, and the 
LB prediction can therefore depend non-trivially (up to 
a rescaling) only on one parameter. A convenient choice 
of this parameter, which we will use here, is 



E-El 
Ed 



(6) 



where Ejj{M, fo) is the energy of the system with mass 
M and phase space density /o in the degenerate limit, 
i.e., is the energy of the system above the degenerate 
limit normalized to the lowest energy possible for the 
same mass and phase space density. The expression for 
Ed is given in Appendix ^ (see also 0). 

We next describe how we derive, once /3 and fi are 
known, the LB predictions for the various quantities we 
will measure in our simulations. 



A. Spatial distribution 



Using (||), the Poisson equation gives, 



d'^(p{x) 
dx 



2gp{x) = 2g 



fo 



dv (7) 



where p{x) is the mass density profile (which we will refer 
to simply as the "density profile"). It is simple numeri- 
cally to solve this (second order differential) equation for 
(p{x), and then to determine the mass density profile, us- 



ing the boundary conditions ^\x=o = and Lp\x=o — 

B. Velocity distribution 

The velocity distribution may be written 

fo 1 



0iv) 



where 



1 _|_ e'3(V+v-A') a{(p) 
dip{x) 



dip 



a{x) 



dx 



(8) 



(9) 



is, up to a sign, the gravitational acceleration. Using the 
Poisson equation we have 



d''^{x) _ Idia'i^)) 



dx^ 



dtp 
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and therefore 



a{p) ^ ]l'^9 J P{v')dtp'. 



(10) 



(11) 



Using the previously determined p{tp) we obtain 6{v) us- 
ing (|). 



C. Energy distribution 

The distribution of particle energies is defined by 



F{e) ^ I S{e-[-- ifiix)]) ■ f{x, v)dxdv (12) 
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with 



J F{e)de = 1. 



Integrating we obtain 



F{e) = D{e)m 



where 
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a{ip) V 



(13) 



(14) 



(15) 



is the density of states at energy e. 

While the results for p{x) and 9{v) do not depend on 
the choice of the zero point of the potential, this latter re- 
sult does. It is straightforward numerically, to use, rather 
than iflx^o 0, 



v^o = ^\x=o = g 



\x\p{x)dx, (16) 



i.e., that corresponding to a pair potential strictly pro- 
portional to the separation between particles. Given that 
a{x) defined in (^ is necessarily positive for all a; 7^ 0, 
this is the minimum value of the potential (and of the 
energy particle energy). Adaptin g th is definition the en- 
ergy distribution is still given by (p^, but with D{e) = 
for e < ipQ and 
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d(f. 



(17) 



D. Order parameters 

In order to characterize and compare the macroscopic 
properties of QSS it is convenient to calculate specific 
moments of the phase space distribution (rather than to 
study always the full distribution). As discussed in |2^ 
a particularly relevant choice can be normalized "crossed 
moments" which give a measure of the "entanglement" 
of the distribution in space and velocity coordinates, by 
considering 



where Ui is the value measured for the particle i. In 
thermal equilibrium the distribution function is separa- 
ble, and so (fia/s = 0. Further it can be shown easily |2^ 
that the thermal equilbrium solution at any energy is the 
unique separable stationary state, i.e., all QSS are non- 
separable. Thus generically we expect these moments to 
be non-zero in a QSS (although any finite number of them 
can in principle vanish without implying separability). 

Here we will use specifically the two moments 4>ii and 
022 to characterize and compare the QSS we obtain in our 
numerical simulations, complemented when necessary by 
examination of the functions derived above and in some 
cases of the full phase space density. Given the LB solu- 
tions determined above (for waterbag initial conditions) 
it is straightforward to calculate numerically the values 
of and (/)22 predicted by LB for this case. These are 
shown in Fig. ^ as a function of the parameter (which, 
as discussed above, can be taken as the single parameter 
on which the LB result depends). We note that both pa- 
rameters are always negative but increase towards zero 
as we go to the non-degenerate limit. Indeed in this limit 
the LB prediction tends to the (separable) thermal equi- 
librium solution. 
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FIG. 1: The "order parameters" 0ii and 022 of the QSS pre- 
dicted by LB theory for waterbag initial conditions, plotted 
as a function of the normalized energy ^d. 
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for non-zero a and /3, where 

^ ^ _ J J uf{x,v)dxdv 
J J f{x, v)dxdv 

estimated in the discrete system with N particles as 

N 



(19) 



(20) 



IV. INITIAL CONDITIONS 

In our numerical study we consider particles dis- 
tributed initially by randomly sampling different classes 
of waterbag initial conditions, i.e., in which the phase 
space density takes the same value /o everywhere it is 
non-zero. Specifically we consider, in order: 

• Single rectangular waterbags (SRW), in which 
the support of the initial phase space density is a 
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FIG. 2: Realizations with iV = 5000 particles of an SRW 
initial condition (left panel) and DRW initial condition (right 
panel). The two configurations have the same value of (up 
to finite corrections). The units used here are specified at 
the beginning of Sec. M below. 




rectangle centred on the origin, i.e., 

fix, v) = foBixo - xMxo + x)e{vo - v)e{vo + v) (21) 

where O is the Heaviside function. As, in the con- 
tinuum limit, the only parameters in the problem 
are then four — /o, xo, vq and the coupling g — 
there is in fact only one relevant parameter charac- 
terizing the system once units are chosen. A nat- 
ural physical choice of this parameter is the initial 
virial ratio i?o, which a simple calculation shows is 
given by 



i?0 



2IL 

Uo 



gMxQ 



(22) 



where Tq and Uq are the initial kinetic and potential 
energies given by 
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Uo = l-gM^xo. 



(23) 



An example of such a configuration with Rq = 0.5 
is given in the left panel of Fig. ||. 

As discussed above the LB prediction also depends 
on only one parameter, which we can take to be ^u, 
the ratio of the energy of the configuration to that 
of the degenerate limit of LB (i.e. the minimum 
allowed energy of the given mass at phase space 
density /o). The energy and mass in the limit of 
a de gener ate s ystem are given as functions of fi by 
(B5) and (B4). Eliminating /i we obtain 



En = 



123 



-gxoM'R§ 



(24) 



and thus 
E 

^E^~~ 



123 



% = 0.688(^ 



Rn 



(25) 

where i?(|, |) is a beta function. This expression 
is is plotted in Fig. |. The SRW with i?o = 1 is 
thus the lowest energy configuration, and there are 
otherwise two values of Rq for each value of £,d- 



FIG. 3: as a function of Ro for a SRW initial condition. 



• Double rectangular waterbag (DRW), in 

which the support of the initial continuum phase 
space density is like that shown in the right panel 
of Fig. I 

f{x,v) = foQ{x + xi)e{xi - x)eiv + vi)eivi - v) 
+ foQix + a;2)e(-a;i - x)Q{v + V2)Q{v2 - v) 

+ fo&{x - X2)&{XI - X)e{v + V2)6(W2 - v) . 

As this has two additional parameters compared to 
the SRW, it is effectively a three parameter family 
of initial conditions, which coincides with the SRW 
when vi — V2, xi = or xi = X2- When they differ 
from the SRW, they are spatially inhomogeneous, 
with a ratio of densities 5 — — vtv the two different 
regions. We will choose to characterize them by 
this parameter, together with and the initial 
virial ratio i?o. LB theory thus predicts that the 
final state should be independent of i?o and 5 at 
given The relevant expressions for the kinetic 
and potential energies of the DRW configuration 
are given in Appendix 

• Disjoint waterbags (DW), in which the initial 
phase space density is made of two disjoint regions 
with simple shapes, either rectangular or elliptical. 
We will use such configurations to further explore 
some of the conclusions draw from the study of the 
SRW and DRW configurations. 



V. NUMERICAL RESULTS 

A. Choice of units 

Unless otherwise indicated our results will be given in 
units fixed by taking g — 1, AI = 1, and Lq — \ where 
Lo the initial linear size of the system, i.e., the distance 
between the outer extremities of the theoretical waterbag 
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initial condition. This implies that the unit of time is 

1 



tc = 



(26) 



where po = M/Lq is the initial mean mass density. This 
is simply a characteristic time scale for a particle to cross 
the system. In the cold limit (i.e. with zero initial ve- 
locity, with Rq — >■ 0) of the SRW initial conditions, it 
corresponds exactly to the time in which all the mass 
falls to the centre of the system. 



B. Attainment of QSS and their characterization: 
generalities 

That the SGS model with a large number of particles 
— just as such 3D self-gravitating and other long range 
interacting systems which have been studied in the litera- 
ture — give rise typically to QSS starting from initial con- 
ditions such as those above has been discussed elsewhere 
in many studies (see references given in the introduction) . 
The attainment of a QSS should be tested, in theory, by 
considering the full phase space density coarse-grained 
at some chosen scale. One would then verify whether its 
evolution after some initial period (of violent relaxation) 
tends to 



fix, V, t) = fqssix, v) + 5f{x, v, t) 



(27) 



where the amplitude of the fluctuations \5f{x,v,t)\ de- 
creases as N increases. For our study here, in which we 
consider how the properties of these QSS depend on the 
initial conditions, what is of importance is that we evolve 
the corresponding system to a time at which the approx- 
imation (|2^) indeed holds well, for N sufficiently large 
so that the fluctuations Sf{x,v,t) introduce a negligible 
uncertainty into the quantities used to characterize the 
QSS. 

In practice numerical limitations on N make a direct 
analysis extremely difficult, and one typically considers 
the behavior of single macroscopic parameters, such as 
the virial ratio, or the magnetization in models (e.g. the 
HMF model) where it is defined. This is then comple- 
mented by a visual inspection of the system represented 
in phase space. To describe the properties of the QSS 
one then considers typically the density profiles, veloc- 
ity and/or energy distribution. We have shown in [ pT| , 
where we studied the very long time behaviour of QSS re- 
sulting from SRW initial conditions, that the parameters 
and (j)22 defined above are very useful macroscopic 
"order parameters" , which can be used to diagnose both 
the attainment of a QSS and to characterize this state. 
We will use them here for the same purpose, supplement- 
ing their calculation where necessary, or interesting, by 
a fuller analysis of the distribution functions. 

To determine whether a QSS is reached, and on what 
time scale, we thus study firstly the evolution of the virial 
ratio, and of 0ii and 022- While the characteristic time 
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FIG. 4: Temporal evolution of virial ratio (top panel), and 
(fill and <j}22 (lower panel) starting from a realization with 
N — 10* particles of a DW initial condition (shown in first 
panel of Fig. ^ . The time units here are such that tc = a/s/S, 
i.e., t = WKi 8.2ta. 



for the mean-field dynamics is of order tc defined above, 
the completion of relaxation to QSS (in the sense de- 
fined above) takes typically of order several tens to one 
hundred tc for SRW initial conditions. Further this time 
depends, unsurprisingly, on the nature of the initial con- 
dition, with very cold initial conditions — further from 
virial equilibrium initially — taking significantly longer 
to relax. 

For DRW and DW initial conditions we observe even 
greater variation in the time for full relaxation to a QSS 
than for SRW, with, in some cases, significant persistent 
fiuctuations in the macroscopic parameters. An exam- 
ple of such a case is shown in Fig. |4| in which the upper 
panel shows the evolution of the virial ratio and the lower 
panel that of the parameters (pu and (j>22, for a DW ini- 
tial condition sampled with = 10"* particles. The full 
phase space plot is shown in Fig. ^. This reveals that it 
is a persistent "rotating hole" feature in the phase space 
which gives rise to the (small but clearly visible ^coherent 
fluctuations in the averaged parameters in Fig. This is 
precisely the kind of effect which has been documented 
in the two studies |2^ mentioned in the introduction, 
and which has been argued in this context to show that 
LB theory is incorrect (as it predicts, by construction, the 
attainment of a time independent phase space density). 
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ever, would be required to establish this conclusion more 
definitively for a broader range of initial conditions, and 
to exclude notably that collisional relaxation may play 
some role. 



C. SRW initial conditions 
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FIG. 6: The density profile (top left), velocity distribution 
(top right) and energy distribution (bottom) for the QSS ob- 
tained starting from SRW with Ro — 0.1. The solid lines are 
the corresponding LB predictions. 
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FIG. 5: Phase space plot of particle trajectories evolved 
from the DW initial condition shown in the first panel (with 
= 10* particles). The time units are those indicated in the 
previous figure. 



While the hole we observe is clearly visible at t = 500 
and indeed rotates in phase space, the subsequent two 
panels show that it slowly disappears on a time scale of 
order a few thousand dynamical times. Thus it appears 
that the relaxation of these holes simply represents a pro- 
longation of the collisionless relaxation to a well defined 
QSS, as no tendency of the system to evolve towards ther- 
mal equilibrium (corresponding to (j^n and (f>22 equal to 
zero) is evidenced on this time scale. Further study, how- 



The density profiles, velocity distributions and en- 
ergy distributions in the QSS obtained starting from 
SRW configurations with Rq — 0.1,0.5,1 are shown in 
Figs. H, and |[ These correspond to averages over 
30 realizations of each initial condition sampled with 
N — 5000 particles, taken at t = 200ic, by which time 
the QSS is well established. In each case the LB pre- 
dictions given in Sec. HI are shown also, correspond- 
ing to — 0.56, 0.08, 0.03 respectively. As observed 
already in early studies |l^, ^ the qualitatively most 
striking deviation from the prediction of LB theory is 
marked by the appearance of a "core-halo" structure, 
most clearly evident in the energy distribution obtained 
from the Rq — 0.1 initial condition. On the other hand, 
as underlined in the more recent study of [ p^ for these 
same initial conditions, the agreement of the LB theory 
with the observed QSS is in fact quite good for the case 
Rq — 1, while the discrepancy progressively increases as 
Rq deviates from unity and a core-halo type structure 
appears. 

Shown in Fig. ^are the values of the parameter (f>ii and 
022 in the QSS, and the values predicted by LB theory. 
This plot summarizes in a simple manner the conclusions 
above: the theory works quite well quantitatively at the 
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FIG. 7: The density profile(top left), velocity distribution 
(top right) and energy distribution (bottom) for the QSS ob- 
tained starting from SRW with Ro — 0.5. The solid lines are 
the corresponding LB predictions. 



We note that these single parameters, 0ii and 022, 
actually allow a better diagnosis of the closeness to LB 
theory than the examination of the full density or veloc- 
ity distribution functions. Indeed comparing just these 
two latter functions with the LB predictions, we might 
conclude that the agreement is almost perfect. The en- 
ergy distribution, on the other hand, allows one to see 
clearly the discrepancies, which are then reflected well in 
(pii and 022 ^ When considering a larger space of initial 
conditions, as we do now, it is very convenient to use 
these parameters as diagnostics of the validity of LB. 
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FIG. 9: 011 (top) and 022(bottom) in QSS as the function of 
Ro- The line indicates 0^/3 calculated by LB stationary state 
and the circle is the value obtained by numerical simulation. 



FIG. 8: The density profile(top left), velocity distribution 
(top right) and energy distribution (bottom) for the QSS ob- 
tained starting from SRW with Rq = 1. The solid lines are 
the corresponding LB predictions. 



lowest energy state corresponding to _Ro = but devi- 
ates greatly as we go towards the less degenerate initial 
states. Further the plot shows that the theory gives very 
qualitatively the correct behavior of the parameters — 
they increase monotonically with the initial ^d- At low 
degeneracy the sign of these parameters is a result of the 
formation of a core which is colder than predicted: there 
is in this case an excess of low velocity particles at small 



D. DRW initial conditions 

As described above the DRW initial conditions allow us 
to test further a basic prediction of LB theory: the same 
QSS should result starting from any initial configuration 



^ This "efficiency" of these parameters as diagnostic tools was 
noted in [ ^ , where it was shown, notably, that they could iden- 
tify clearly stationary states arising from certain initial condi- 
tions as QSS rather than the thermal equilibrium states which 
previous studies tea] had mistakenly inferred them to be based 
on an analysis using p{x) and 0(v). 
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in the range of accessible "microstates" at given mass and 
energy. For ID gravity and waterbag initial conditions, 
this means the QSS obtained should be the same at a 
given independently of the shape of the waterbag. As 
discussed the DRW gives us a two dimensional space of 
such configurations, which we choose to parametrize by 
the initial virial ratio Rq and density contrast S. 

For each of the three values of corresponding to the 
SRW initial conditions above, we have simulated twenty 
different initial conditions chosen to explore the available 
{Rq, S) space. In each of Figs. |lO|, |l^ and ^ are shown 
two plots: one shows the initial conditions in the {Rq, S) 
plane at the given value oi£,D, the other the QSS obtained 
from them as represented in the plane {(fiii, 4>22)- The 
results are, as for the SRW above, averages over 30 real- 
izations of each initial condition sampled with N = 5000 
particles, taken at t w 200tc. The fact that the spread in 
values of (i?o, 5) is much smaller at smaller is simply 
a reflection of the fact that as one goes towards the de- 
generate limit = the constraints limit the possible 
deformations more and more. 
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FIG. 11: Same as in Fig. |l^, but for DRW initial conditions 
with = 0.08, i.e., equal to that of SRW with Ro = 0.5. We 
keep the scale as in the previous figure for easier comparison. 
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FIG. 10: The upper plot represents the twenty different DRW 
initial conditions with = 0.56 (i.e. equal to that of SRW 
with Ro = 0.1) according to their values of Ro and 5. The 
lower plot represents the values of (022, 0ii) measured in the 
resulting QSS. The LB prediction lies at the centre of the 
small circle. The unfilled points in the upper plot correspond 
to the four initial conditions which give QSS closest to the LB 
prediction. 

In continuity with what we observed for the SRW, the 
results show that LB theory works reasonably well at 



the two lower values of — the QSS varies only very 
little over the range of different initial conditions — it 
is grossly violated as we go towards the non-degenerate 
limit. Indeed the order parameters for QSS obtained 
starting from the same S^d can differ in sign. Direct 
analysis of the distribution functions confirms that this 
corresponds to QSS which are completely different. On 
the other hand, certain initial conditions at = 0.56 
— those corresponding to the unfilled points in the up- 
per plot of Fig. ^ — do appear to give QSS close to 
the LB prediction. To assess whether this is really the 
case the density profiles, velocity and energy distribution 
functions for two of these are shown in Figs. |l^ and 
While the agreement with the theoretical curves is not 
perfect, it is comparable than that obtained for the ini- 
tial conditions with = 0.03 — indeed the discrepancy 
between the LB prediction and the observed distribu- 
tions is no more than observed above for the SRW initial 
conditions with — 0.03. 

The strong deviations from the LB prediction, just as 
in the SRW, manifest themselves in the shift towards pos- 
itive values of <j)ii and (/)22 • Direct inspection of the distri- 
bution function of energy shows that this reflects again in 
all cases the appearance of a pronounced core-halo type 
structure. Inspection of the plot of the initial conditions 
in the (i?o, 5) space for S^d = 0.56 shows that all the cases 
which approach LB (unfilled points) are characterized by 
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FIG. 12: Same as in Fig. |l^, but for DRW initial conditions 
with = 0.03, i.e., equal to tliat of SRW with iio = 1- Same 
scale as in Fig. hd for ease of comparison. 



an initial virial ratio near unity, while the density con- 
trast parameter S appears to be irrelevant. On the other 
hand i?o ~ 1 is clearly not a sufficient condition to guar- 
antee agreement with LB. 

These results suggest therefore that LB theory works 
reasonably well always near the degenerate limit, and 
also for much higher energies for very specific initial con- 
ditions. In these cases, which seem to correlate strongly 
with an initial virial ratio near unity, the formation of 
a core-halo structure, not predicted by LB theory, is 
avoided. 



E. DW initial conditions 

To further explore these findings, and in particular to 
investigate the relevance of the initial virial ratio as a 
parameter, we consider finally a few other "disjoint" wa- 
terbag initial conditions as described above. We report 
results for the four cases shown in Fig. ^ Each of the ini- 
tial conditions has been adjusted to have i?o = 1, and the 
values of the normalized energy are = 1.59,0.58,0.49 
and 0.23 for DWl to DW4 respectively. We take in each 
case a single realization with N = 10^ particles, and cal- 
culate a time average by sampling on 100 equally spaced 
time slices in the time window [4000, 5000] (in the time 
units of our simulation, which differ in each case from 



1.2 - 




-1.6 -0.8 0.8 1.6 -1.4 -0.7 0.7 1.4 




0.5 1 1.5 



e 

FIG. 13: Density profile(top-left), velocity (top-right) and en- 
ergy distribution (bottom) for DRW initial conditions with 
= 0.56 (i.e. the same energy as the SRW with Ro — 0.1), 
iio = 1.39 and S = 0.054. 
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FIG. 14: Density profile(top-left), velocity (top-right) and en- 
ergy distribution (bottom) for DRW initial conditions with 
= 0.56 (i.e. the same energy as the SRW with Ro — 0.1), 
Ro = 1.017 and 5 = 9.861. 



units with tc = 1 by a numerical factor of order unity) . 
Shown in Fig. |l6| are the QSS obtained as represented 
in the ((/)ii,022) plane. In each case the filled symbol 
represents the corresponding LB predictions. Compar- 
ing to the results for SRW and DRW initial conditions, 
the QSS appear in all cases much closer to the LB predic- 
tions. This is confirmed by inspection of the distribution 
functions, which are shown for DW2 in Fig. |l^ and for 



11 



1.5 
1 

0.5 


-0.5 
-1 
-1.5 



1 


c 


iA/1 


1 







1.5 
1 

0.5 


-0.5 
-1 
-1.5 





DW2 









-0.5 



0.5 



-0.5 



0.5 



-1.5 -1 -0.5 0.5 1 1.5 



-1.5 -1 -0.5 0.5 1 1.5 



1.5 
1 

0.5 



-0.5 
-1 

-1.5 



1 


DW3 










1 


1 



-0.5 



0.5 




FIG. 15: Four disjoint waterbag initial conditions with the 
number of case indicated in the panel. The corresponding 
are equal to 1.59,0.58,0.49 and 0.23 for case 1 to 4 respec- 
tively. 




FIG. 17: The density profile (top left), velocity distribution 
(top right) and energy distribution (bottom) for the QSS ob- 
tained starting from the DW2 initial conditions (^_d ~ 0.58). 
The solid curves lines are the LB predictions. 



DRl in Fig. |l8|. For the former case the results are as 
close to the LB predictions as for the SRW and DRW 
cases which gave best agreement with LB, with the small 
deviation being visible again in the energy distribution 
but very difficult to discern in p{x) or 9{v). The results 
for the cases DR3 and DR4 are similar. For DRl, on the 
other hand, the deviation from LB is much more marked, 
and we see in the energy space that this deviation is asso- 
ciated to the formation of a (in this case very small) core. 
Very much in line with the results for SRW and DRW ini- 
tial conditions, the agreement with LB thus deteriorates 
as one goes away from the degenerate limit. 
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FIG. 16: 022 and (pn of the QSS obtained from the initial 
conditions in the previous figure. The unfilled symbols cor- 
responds to the values obtained from numerical simulations, 
while the filled symbols are the LB predictions. 



In summary these results confirm the conclusion drawn 
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FIG. 18: The density profile (top left), velocity distribution 
(top right) and energy distribution (bottom) for the QSS ob- 
tained starting from the DWl initial conditions (^_d ~ 1.59). 
The solid curves lines are the LB predictions. 



from the analysis of the SRW and DRW waterbags: the 
LB predictions are always reasonably good — and excel- 
lent for the spatial and velocity distributions - for (wa- 
terbag) initial conditions with low , but even at higher 
values good agreement can be obtained in cases charac- 
terized by an initial virial ratio of order unity. Further 
deviation from LB is always characterized by the appear- 
ance of a core-halo type structure. 
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VI. COMPARISON WITH THEORETICAL 
PROPOSALS BEYOND LB: DIRECT ANALYSIS 
OF PHASE SPACE DENSITY 

Let us consider how well two recent proposals in the 
literature can account for the properties of the QSS we 
observe: 

• Yamaguchi ^ studies the SGS model for SRW 
initial conditions, and notes (as was remarked also 
in early studies 
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that the breakdown of LB 
theory is associated with the apperance of a core- 
halo structure. He proposes a phenomenological 
adaptation of LB theory which he uses to fit the 
resultant core, in which the LB theory is applied 
only to the mass and energy associated to the core. 
In practice this means that one parameter is mea- 
sured a posteriori from the observed QSS. 

• Levin et al. in a series of works on other mod- 
els — plasmas Q, 3D gravity [Q, 2D gravity Q 
and, most recently, the HMF model ||ll[ — have 
proposed that, when LB theory breaks down, QSS 
correspond to the phase space density: 

fix, v) = foMeF - e) + X e(e - 6^)9(6,, - e)] . (28) 

As in the case of |l^ , this involves the addition of 
one parameter compared to LB theory. However, a 
physical explanation is proposed for the core-halo 
form of (^8|), and a prediction for this additional pa- 
rameter is derived from the initial conditions: An 
analysis of particle dynamics in the coherent os- 
cillating field associated with the relaxation shows 
that there are dynamical resonances which allow 
particles to gain energy, with Ch corresponding to 
the maximal energy which can be attained in this 
way. Assuming that resonance effect is "shut off' 
only by the upper bound on the phase space density 
imposed by the coUisionless dynamics, the ansatz 
( p8| ) is the simplest one possible for the QSS which 
will result. 

To evaluate the validity of these approaches in this 
model, we consider directly the measured phase space 
density, /(e), as a function of particle energy. To do so we 
measure the (averaged) values of the potential (l){x) and 
a((/?) in the QSS, and then use ( p^ ) to calculate the phase 
space density D{e). Shown in Fig. ^ are the results for 
eight chosen cases from the DR W in itial conditions with 
= 0.56 considered in Sec. VD above. In Fig. ^ a 
plot of exactly the same data is given, but now displaying 
the logarithm of the absolute value of //(/o — /) as a 
function of e (which in LB theory gives a straight line 
with slope (3). Our choice is representative of the whole 
batch of initial conditions, in that 1) most QSS have a 
clear core-halo structure, and 2) those that do not agree 
reasonably well with the LB prediction. Indeed the two 
configurations in the uppermost panel of Fig. n9 are the 



same two cases for which the full distribution functions 
were shown in Figs. ^ and |lj. 

In Fig. |l^ a vertical line indicates the initial phase 
phase density /o, so that it is clear that whenever a 
core appears it is indeed degenerate. While the mea- 
sured phase space distributions are clearly more struc- 
tured than (|2^), in most cases this simple ansatz gives 
a reasonably good fit (i.e. about as close to the phase 
space density as the LB profile is to the observed one 
in the cases where it has been considered to work well 
above). The slightly greater structuration of the phase 
space density compared to the ansatz of (p^ ) can also be 
seen in Fig. ^0|, which shows in particular that the dif- 
fuse halo, when present, although close to flat, appears 
clearly more consistent with a Maxwell-Boltzmann form 
(i.e. the non-degenerate limit of LB theory). In this re- 
spect we note that Levin et al. have not tested their 
ansatz directly against the phase space density, but have 
used it to derive predictions for p{x) and 0{v) which have 
been compared with those observed. As we have seen in 
comparing numerical results with LB predictions above, 
these quantities typically wash out structure in energy 
space and make it difficult to see discrepancies which are 
localized in this space. We note further that our finding 
that it is initial conditions with i?o ~ 1 which suppress 
the core-halo formation — and lead to QSS in reason- 
able agreement with LB — appears completely coherent 
with the mechanism described by Levin et al.: when the 
system starts close to virial equilibrium, the relaxation 
is typically indeed much "gentler", simply because the 
system does not undergo the large contractions and ex- 
pansions which result necessarily if there is a large im- 
balance between the initial potential and kinetic energy. 
It is precisely such macroscopic oscillations of the system 
which drive the resonances analyzed by Levin et al.. 

We consider finally comparison of our results with an 
analytical treatment of coUisionless relaxation developed 
in |2(| (see also further references therein). This work 
develops, under certain approximations and hypotheses, 
a kinetic equation for coUisionless relaxation — similar to 
the Lenard-Balescu equation for collisional relaxation — 
with a term describing relaxation towards the LB equi- 
librium. One feature of this term is that it involves an 
effective space and time dependent diffusion coefficient, 
which is proportional to the product /(/o — /). Thus the 
theory suggests that relaxation should be expected to be 
most inefficient when / is close to degenerate (/ ~ /o) 
or very small (/ ~ 0). In regions of energy where relax- 
ation is more complete, the distribution is expected to ap- 
proach the LB form, but with values of the parameters /3 
and fi different from those in the global LB equilbrium. 
Our results in Figs. |l^ and ^ do appear to be quite 
consistent with these qualitative predictions: indeed this 
theory would appear to account for why it is core-halo 
type states, whose dynamical origin is explained by Levin 
et al., which do not relax to the (global) LB equilibrium. 
In all cases the results in Fig. show a region where 
the halo distribution is very consistent with a Maxwell- 
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FIG. 19: Phase-space density as a function of energy /(e) for 
eight representative cases of DRW initial conditions with = 
0.56 (corresponding to Ro — 0.1 for SRW). The two upper 
panels correspond to the two cases for which the distribution 
functions are shown in Figs. |l|and|l^ where the QSS is close 
to LB. The dashed horizontal line indicates the initial phase 
space density, fo, while the continuous lines correspond to the 
LB prediction. 



FIG. 20: Exactly the same data as in the previous figure, but 
now with the logarithm of the absolute value of f /{fo — f) 
plotted as a function of particle energy e. The dashed lines 
represent the predictions of LB theory (which become straight 
lines of slope /3 in this representation). 



VII. CONCLUSION AND DISCUSSION 



Boltzmann form with an inverse temperature lower than 
that of the global LB prediction (dashed line). This can 
be interpreted, following , as a "mixing region" where 
the (in this case, non-degenerate) LB distribution applies 
locally, while the deviation from the (local) LB form at 
higher and lower energies is considered as due to incom- 
pleteness of relaxation in these regions. Further the fact 
that the observed distributions are, compared to the ex- 
trapolated straight line ("local" LB) fit in the "mixing 
region" , sensibly higher at lower energies and lower at 
the highest energies is also in apparent agreement with 
the kinetic theory described in . 



We summarize now our principal conclusions: 

• Attainment of QSS in the SGS model: In all 

cases we have considered QSS do appear to be at- 
tained, but the time-scales for relaxation to them 
can vary very considerably. In some cases rotat- 
ing "holes" formed in the phase space density dur- 
ing the initial phase of violent relaxation (^ lO^tc) 
survive for quite a long time, disappearing only on 
times scales of order (~ lO'^ic). As on these latter 
time scales the system shows no apparent tendency 
to relax towards its thermodynamic equilibrium, we 
conclude that this is simply a manifestation of slow 
collisional relaxation, and does not imply that QSS 
are not attained as argued by A fuller 
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study of the possible N dependence of such relax- 
ation would be useful to establish this conclusion 
more firmly (but would be numerically challeng- 
ing). 

• LB theory in the SGS model: As was clear 
already from early studies and confirmed 
by more recent ones such as [p_5[, LB theory is not 
an adequate theory for understanding fully, or even 
approximately, the properties of QSS arising from 
violent relaxation in the SGS model for arbitrary 
initial conditions. However, it is by no means an 
irrelevant theory to understanding these QSS. Our 
study of a quite broad range of initial conditions 
shows that the space of QSS in this model divides 
quite neatly into two: those for which LB works to 
quite a good approximation, and those for which 
the phase space density is characterized by a de- 
generate core, taking a form generally quite close 
to the simple ansatz ( p8|) proposed by Levin et al. 
The initial conditions in the former class are either 
close to the degenerate limit, or in other cases char- 
acterized by initial virial ratio of order unity. These 
conditions are precisely those, in line with what 
has been described by Levin et al., which suppress 
resonances which otherwise act very efficiently to 
produce the degenerate core-halo structure. 

• Accuracy of predictions of QSS in the SGS 
model: While, as just described, the QSS which 
result from violent relaxation divide into those 
which are close to the LB theory, on the one hand, 
or to the ansatz of Levin et al., on the other, the 
accuracy of the associated predictions is at best ap- 
proximate: in no case do we see a perfect agree- 
ment with either LB theory or the ansatz (p^). We 
underline that in this respect the spatial distribu- 
tion of mass p{x) and velocity distribution 9{v) are 
rather poor tools for diagnosing the agreement be- 
tween observations and theory, as they wash out 
deviations which are most pronounced in energy 
space. We have also noted the apparent coherence 
of our results with the qualitative predictions of the 
kinetic theory approach described in p6| . 

Numerous results in the literature on various other 
models (see references in introduction) for specific ranges 
of initial conditions suggest that these latter two conclu- 
sions, and probably the first also, might apply much more 
generally to long-range systems. Further detailed investi- 
gation of such models, and in particular of broader classes 
of waterbag initial conditions like those considered here, 
or, for example, "multi-level" waterbag initial conditions 
would be required to establish if this is the case. 

For the SGS it would be interesting to apply the analy- 
sis described by Levin et al. to determine a prediction of 
the form ( p8| ) for different initial conditions, and see how 
well it does in approaching the observed QSS. In this 
respect it is interesting perhaps to note that, at given 



value of this is a one parameter family of solutions, 
so that it predicts QSS lying on a curve in the (^iii, 4>22) 
plane. In Fig. |o[ we see that the QSS obtained from 
the two parameter family of initial conditions at a fixed 
£_D = 0.56 do approximately collapse onto a curve. We 
would expect the degree to which the simple ansatz ( |2^ ) 
can fit the QSS to be well characterized by determining 
the prediction it gives in this plane. 

Of particular interest is of course the original context of 
3D self-gravitating systems, to which the initial study of 

fl for SRW suggests these conclusions may indeed apply, 
s mentioned, however, the results reported have been 
based, in this case, on examination of the density profile 
p{x) alone, while the energy distribution is probably a 
finer diagnostic tool as we have seen here. In forthcom- 
ing work we will study this case, and discuss the possible 
relevance of our findings in the astrophysical context. In 
this respect we note one of the reasons why LB theory has 
not played — at least for what concerns it detailed pre- 
dictions — a role in astrophysics is that these predictions 
depend on unobservable initial phase space densities. In 
contrast the prediction of a degenerate core in many cases 
would give a simple link between observations and initial 
conditions, which may be of practical relevance notably 
in constraining the parameters in theories of structure 
formation in the universe. 

The simulations were carried out in large part at the 
Centre de Calcul of the Institut de Physique Nucleaire et 
Physique des Particules. We are particularly grateful to 
Laurent Le Guillou for advice and help on use of these 
computing resources. We thank B. Marcos for useful dis- 
cussions, and P.H. Chavanis for many helpful remarks, 
and in particular for suggesting the plot in Fig. 



Appendix A: Determination of /3 and fj, 

In general /3 and fi cannot be calculated analytically, 
so we solve for them numerically as follows. The mass 
normalization (Q) condition is 

/OO POO 
/ f{x,v)dxdv. 

Integrating over v, and changing the coordinate x to (p{x) 
just as in (|^), we obtain 

M=Ar r i^dvd^. (Al) 

Jo Jo a{ip) 
The total energy constraint (||), i.e., 

roc i-oc, 2 / \ _ 

^ = / / {^ + ^)IMdxdv, 

J —OO J —OO ^ ^ 
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can likewise be rewritten as 
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oo poo 2 





oo pOC 



a{ip) 



dvdip 
dvdip 



(A2) 



= T + C/, 



where T is total kinetic energy and U is total potential 
energy. We can then use the virialization condition, 2T = 
U, to obtain 



E = 6 



OO poo 2 



a{(p) 



dvdip. 



(A3) 



The determination of the parameters /3 and fi in the LB 
solution (^) can then be cast as the problem of finding 
the solutions of the equations 

F{p,fi) - 
G(/3,A*) = 0. 



where 
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(A4) 
dvdip. 
(A5) 



Following a standard method we write the matrix equa- 
tion 
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(A6) 



where dF and dG denote the infinitesimal changes of F 
and G when (/3, /i) change to {fH + dp, fi + dfi), we start by 
guessing a pair of (/?, fi) and then determining the new 
{P',fi') = {p + Al3,fi + Afi) using 
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dfi 


dti 


dG 


dG 


dp 


dfi 



AF{(3,fi) 
AG{p,^l) 



(A7) 



where {AF,AG) = {-F{(3, n), ~G{l3, fi)). We then iter- 
ate until AF and AG converge to 0. With a reasonable 
guess for the starting values of /3 and /i, good conver- 
gence is attained within a few iterations, as illustrated in 
Fig. El] for a typical case. 



Appendix B: Degenerate limit of LB theory 

For completeness we reproduce here the analytic re- 
sults of 1l4| (see also [E9|) for the degenerate limit of 
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FIG. 21: The values of /? and /j, obtained in the successive 
steps of our iterative numerical calculation, for a typical case. 
The units here are those used in our numerical calculation 
{M — N , g = 1 and Lq — TV), different to those in which our 
results in the main text are given. 



the LB distribution function 
/? — > oo, in which 



fix,v) 

The density profile is then 



fo, e{x,v) < n 
0, otherwise. 



This corresponds to 
(Bl) 



p{ip) = 2 / fodv 
Jo 

= 2V2Mfi-p]i 
and therefore, using (O), 



aiip) = A-2i-{^)^^^-i^-p)i]i. 

The mass normalization then yields 
Jo 

4^2/0 /-^ {fi-^)-2d^ 



(B2) 



(B3) 



^/ 16V|g/o Jq {fl2 - {p. - p)2)2 

which can be integrated to give 

M = 2f(f)iMi 



(B4) 



Using the expression ( AS ) for the total energy we have 

fM /-V^Im-v) ,,2 



6/o 



a{ip) 



dvdip. 



16 



Integration first over v gives 



netic energy is 



Ed - 24(^)2 / 3 , ,3,1 rfv^ 
and then, on integrating by parts, we obtain 



To = M 



2{x2 — xi)v^ 2xivl 



3 35 Jo 

where ip' — i^^-^)^ ■ The integral can be expressed as a 
beta function, and the result can thus be written 



3 ^3g^ 4' 3^^ 



(B5) 



It is simple to show from (A2) that 4^ < in general, 



tending asymptotically to as /3 — ?> cxi. This is thus in- 
deed the minimal possible energy corresponding to given 
M and fo- 



Appendix C: Generation of DRW initial conditions 



For the DRW phase space density defined in Sec. EV 
a direct calculation gives immediately that the initial ki- 



and the initial potential energy 



Ax 



2x1 



+8foViV2g[xi{xl - xl)] + ^fovfgxl 



These indeed reduce to the corresponding expressions 
(p3| ) for the SRW (when we set Xi = 0, xi = X2 or 
vi — V2)- To g enera te the specific initial condition re- 
ported in Sec. VD^ we do a random sampling in in 
xi,X2, vi and V2 at fixed fo , E and M fixed (which 
implies that is fixed). We then choose configurations 
with Ro = 2To/Uo and S as various as possible. 
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